///\todo upgrade

//////////////////////////////////////////////////////////////////////////////////
//                                                                              //
//  This file is part of the buola project (https://code.google.com/p/buola/).  //
//                                                                              //
//  Copyright(c) 2007-2012 Xavi Gratal                                          //
//  gratal AT gmail DOT com                                                     //
//                                                                              //
//  Buola is free software: you can redistribute it and/or modify               //
//  it under the terms of the GNU General Public License as published by        //
//  the Free Software Foundation, either version 3 of the License, or           //
//  (at your option) any later version.                                         //
//                                                                              //
//  Buola is distributed in the hope that it will be useful,                    //
//  but WITHOUT ANY WARRANTY; without even the implied warranty of              //
//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               //
//  GNU General Public License for more details.                                //
//                                                                              //
//  You should have received a copy of the GNU General Public License           //
//  along with buola.  If not, see <http://www.gnu.org/licenses/>.              //
//                                                                              //
//////////////////////////////////////////////////////////////////////////////////

#ifndef _BUOLA_MAT_INITIALIZERS_H_
#define _BUOLA_MAT_INITIALIZERS_H_

#include <buola/mat/detail/generators.h>
#include <buola/mat/detail/cuniforminitializer.h>
#include <buola/mat/detail/clinearinitializer.h>
#include <buola/mat/detail/cidentityinitializer.h>

namespace buola { namespace mat {

inline detail::CUniformInitializer<detail::COnesGenerator> ones(size_t pRows,size_t pCols=1)
{
    return {pRows,pCols,{}};
}

inline detail::CUniformInitializer<detail::CZerosGenerator> zeros(size_t pRows,size_t pCols=1)
{
    return {pRows,pCols,{}};
}

template<typename tConstant>
inline detail::CUniformInitializer<detail::CConstantGenerator<tConstant>,tConstant> constant(size_t pRows,size_t pCols,
                                                                           tConstant pConstant)
{
    return {pRows,pCols,{pConstant}};
}

template<typename tDistribution=std::uniform_real_distribution<double>>
inline detail::CUniformInitializer<detail::CRandomGenerator<tDistribution>>
random(size_t pRows,size_t pCols,const tDistribution &pDist=std::uniform_real_distribution<double>(0,1))
{
    return {pRows,pCols,{pDist}};
}

template<typename tDistribution=std::uniform_real_distribution<double>>
inline OEnableIfNot<QIsNumeric<tDistribution>,detail::CUniformInitializer<detail::CRandomGenerator<tDistribution>>>
random(size_t pRows,const tDistribution &pDist=std::uniform_real_distribution<double>(0,1))
{
    return {pRows,1,{pDist}};
}

inline detail::CIdentityInitializer<> eye(size_t pRows,size_t pCols)
{
    return {pRows,pCols};
}

inline detail::CIdentityInitializer<> eye(size_t pSize)
{
    return {pSize,pSize};
}

template<typename tScalar>
inline detail::CLinearInitializer<tScalar> linspace(size_t pRows,tScalar pMin,tScalar pMax)
{
    return {pRows,pMin,(pRows==1)?0:((pMax-pMin)/(pRows-1))};
}

template<typename tScalar>
inline detail::CLinearInitializer<tScalar> range(tScalar pMin,tScalar pMax,tScalar pStep=tScalar(1))
{
    return {size_t((pMax-pMin)/pStep+1),pMin,pStep};
}

template<typename tExpr>
int sample(const CExpr<tExpr> &pE)
{
    const auto &lE=detail::element_accessible<detail::ACCESS_FAST>(pE.Derived());
    
    double lCQX=max(lE);
    size_t lN=lE.size();
    
    while(true)
    {
        int lX=random::get<int>(lN-1);
        double lU=random::get<double>();
        if(lU<lE(lX)/lCQX) //only accept it if u<p_x/(c*q_x)
            return lX;
    }
}

/*namespace mat*/ } /*namespace buola*/ }

#endif
